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VARIATIONAL  PRINCIPLES  AND  GENERALIZED  PRINCIPLE  FOR  THREE- 
DIMENSIONAL  TRANSONIC  FLOW  WITH  SHOCK  WAVES  IN  A  ROTATING 
TURBO-IMPELLER 

by  Liu  Gaol i an 

(Shanghai  Institute  of  Mechanical  Engineering) 

Abstract 

In  this  article  a  family  of  variational  principles  for 
three-dimensional  transonic  steady  relative  flow  with  em¬ 
bedded  shocks  in  a  turbomachine  impeller  of  axial  radial 
or  mixed  flow  type  is  developed.  Its  special  feature  is  to 
take  full  advantage  of  natural  boundary  conditions  and 
"artificial  interfaces."  and  besides  converting  all  boundary 
conditions  into  natural  ones,  it  is  also  shown  that  by  taking 
variations  of  the  position  of  unknown  flow  discontinuities 
(such  as  shocks  and  free  trailing  vortex  sheets) ,  all  match¬ 
ing  conditions  across  these  discontinuities,  including  the 
well-known  Rankine-Hugoniot  shock  relations,  can  be  de¬ 
rived  from  the  variational  principles  as  natural  interface 
conditions. 

This  article  is  primarily  intended  to  provide,  in  con¬ 
junction  with  the  discontinuous  finite  element,  a  theoretical 
basis  for  developing  a  new  computational  method  which  allows 
all  flow  discontinuities  (including  the  shock  surface  and 
free  trailing  vortex  sheets)  to  be  singled  out  automatically 
and  distinctly.  Owing  to  the  use  of  a  potential  function,  the 
results  in  this  article  are  limited  to  cases  where  the  Mach 
numbers  before  the  shock  do  not  differ  greatly  from  one. 
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In  the  appendix,  an  analogy  between  the  functional  varia¬ 
tions  with  variable  domain  and  the  Reynold's  transport 
theorem  is  pointed  out  and  used  as  a  basis  for  deriving  a 
general  formula  for  such  variations. 

1 .  Preface 

Following  the  continuous  development  of  high  voltage, 
high  speed,  high  temperature  and  high  efficiency  impellers, 
transonic  flows  increasingly  appeared.  Owing  to  the  simultan¬ 
eous  existence  of  the  subsonic  and  supersonic  areas,  the  flow 
field's  physical  properties  (domain  of  dependence  and  domain 
of  influence)  and  mathematical  properties  (elliptic  type  and 
hyperbolic)  are  all  very  different,  when  the  shape  of  its 
interface  is  unknown,  especially  the  discontinuities  which 
have  yet  to  appear,  the  problem  of  deriving  the  transonic 
flow  field  is  difficult  yet  is  of  real  significance. 

At  present,  the  two  major  methods  for  deriving  the  tran¬ 
sonic  flow  field  are  the  non-steady  variation  method  [3]  and 
relaxation  method.  Although  both  methods  are  automatic,  yet 
neither  can  distinctly  compute  (widened  to  a  difference  net¬ 
work)  the  discontinuities.  Moreover,  the  computation  time  for 
the  non-steady  variation  method  is  very  long  and  it  is  also 
difficult  to  apply.  On  the  other  hand,  although  the  shock 
fitting  method  is  able  to  be  distinct,  yet  it  is  still  unable 
to  compute  the  discontinuities  automatically.  Therefore,  it  is 
very  difficult  to  apply  to  the  problem  of  the  complex  internal 
flow  of  the  shock  system.  The  variational  principles  can  also 
be  taken  as  the  basis  for  deriving  the  transonic  flow  (for 
example,  references  [9-11]),  yet  all  are  limited  by  the 
planar  wing,  especially  the  problem  of  not  taking  into  consid¬ 
eration  the  handling  of  various  unknown  discontinuities. 
Because  reference  [10]  takes  the  positions  of  the  shock  as 
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known,  theoretically  it  is  incomplete  and  unsuitable  for 
application. 

In  the  last  several  years,  the  three-dimensional  aero¬ 
dynamic  theory  for  impellers  has  advanced  rapidly  in  the  area 
of  positive  propositions.  Aside  from  using  the  quasi-three- 
dimensional  theory  which  is  based  on  the  two  relative  flow 
theories  mentioned  in  reference  [1] ,  there  has  been  great 
effort  expended  both  domestically  and  abroad  to  find  a  general 
theory  and  method  for  deriving  three-dimensional  flow  (for 
example,  references  [2-4]). 

This  article,  on  the  basis  of  reference  [14]  and  in  con¬ 
junction  with  the  special  characteristics  of  transonic  flow, 
applies  the  commonly  used  potential  flow  approximation  assump¬ 
tion  and  sets  up  a  family  of  variational  principles  for  three- 
dimensional  transonic  steady  relative  flow  with  imbedded 
shocks  in  an  impeller.  After  combining  this  with  the  discon¬ 
tinuous  finite  element  [5] ,  we  seek  a  theoretical  basis  for 
developing  a  new  computational  method  which  allows  the  various 
unknown  flow  discontinuities  (including  shock  and  free  trail¬ 
ing  vortex  sheets)  to  be  singled  out  automatically  and  dis¬ 
tinctly. 

II.  The  Aerodynamic  System  of  Equations  and  Boundary  Conditions 

of  a  Steady  Relative  Flow  in  a  Rotating  Impeller 

Here  we  will  consider  the  ideal  gas  for  making  a  steady 
three-dimensional  adiabatic  flow  in  a  system  with  a  fixed 
angle  of  velocity  W  rotating  around  a  Z  axis  (see  fig.l). 
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Fig.  1  .  Three  Dimensional  Relative  Flow  Field  in  a  Rotating 
Impeller 

Key:  1.  Entrance 
2 .  Exit 

Due  to  the  fact  that  the  entropy  gradient  after  the  curve 
shock  can  usually  be  considered  very  small  in  the  transonic 
area,  we  can  approximate  that  the  flow  field  has  equal 
entropy  (i.e.  Vs=0).  If  the  gas  entering  in  front  of  the 
rotating  impeller  has  absolutely  no  rotation  (  ^xc^O)  and 

generalized  enthalpy  (or  stagnation  rotor  enthalpy) 
w  2  2 

R=i+  2-w  r  is  homogeneous  (i.e.  vR.=0),then  from  the  gas  flow 
2  2  1 

equation  [1],  wx(  Vxc)  =  VR-T  Vs,  and  the  Thomson  theorem 
we  know  that  the  flow  in  this  rotating  impeller  will  maintain 
absolute  non-rotation  everywhere,  that  is,  VxcSO.  We  know 
from  this  that  there  must  exist  a  velocity  potential  function 
$  (r,  V  ,z)  . 
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1 .  The  Aerodynamic  System  of  Equations 

The  aerodynamic  system  of  equations  for  the  above  mentioned 
flow  can  be  written,  in  the  cylindrical  coordinate  system 
(r,  ,z)  fixed  on  a  rotor,  into  dimensionless  forms  as 

follows  [1] : 


dr 

d<P 

ds 


r  •  V  •  (3A)  -  d&rJ^  +  MMti  +  9(P,Al) 
dr  dip  9m 

—  Ar  -g  -  ^  r  -  (A,  +  A»)r 

ao  &  o  dcp  ft% 

—  -  —  -  A. 
a#  a* 

^  -  ?/*  -  j*.  -  |i  -  J~  (a*  -  a;)}‘ 


o 


(0 

(2) 

(3) 


Here,  besides  the  generalized  enthalpy  modified  use  of  the 
symbol  R,  we  used  dimensionless  quantities  and  symbols  similar 
to  those  of  reference  [8].  We  could  also  use  equation  (2)  to 
eliminate  the  A  in  equations  (1)  and  (3)  and  obtain  two 
equations  which  only  contain  $  and  p  (or  *a>  as  shown  below 
[1]  : 
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In  the  equations 
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Aside  from  this,  we  could  also  substitute  equations  (3*)  and 
(3W)  into  equation  (l1)  and  obtain  an  equation  which  only 
contains  4  . 

2.  Conditions  on  the  Boundaries  and  on  the  Discontinuities 
We  will  only  discuss  the  boundary  conditions  for  whole 
impellers  which  have  not  yet  reached  a  state  of  complex 
blockage . 

Letting  A^  represent  the  section  (area)  of  the  impeller's 
entrance  and  exit  as  well  as  the  entire  solid  wall  surface 
(above  it  there  is  possibly  gas  suction  or  blowing)  of  a  flow 
path  in  an  impeller  (fig.  1) ;  using  Ai  and  A"  to  indicate 

A  ,  \i  Z  /  U 

a  cascade’s  upstream  cycle  reference  boundary  surface,  the 
cycle  distance  between  the  two  is  an  angular  segment  distance 
A*P  -2  If  /N  and  here  N  is  the  blade  number;  using  A^  ^  and 

A^  d  to  indicate  the  cascade's  downstream  cycle  boundary 

surface  (free  trailing  vortex  sheet) ,  the  cyclic  distance 
between  the  two  is  AV  .  Using  A  to  indicate  the  shock, 

S 

S  to  indicate  the  artificial  interface,  and  A  to  indicate 
the  whole  boundary  surface:  A-k^U  a2  UAg(JS  .  The  subscript 
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pr  is  used  to  indicate  the  given  quantity.  Therefore,  the 
boundary  conditions  can  be  expressed  as  follows: 


1)  On  the  A^  surface.  The  given  normal  direction  dense 
flow  distribution. 


HA.  -  <>.)„ 


(4) 


It  shows  that  there  is  a  fixed  distribution  of  blowing  (suction) 
gas  on  the  outside  (wall)  surface. 

2)  On  the  surface.  On  the  upstream  A2  u#  there  is 
cyclicity  of  all  of  the  air  flow  parameter  along  the  perimeter 
and  the  angular  segment  distance  is  the  cycle.  The  vel¬ 

ocity  potential  $  is  then  indicated  with  the  following 
equation: 


<P(r,  <p±\  •  A<p,  *)  —  <P(r,  <p,  *)±*  •  A<P. 

(*-0,1,2,-..)  (5) 

In  this  equation,  -  &  1  is  the  given  value.  On 

U  a  /  U  **  Z  /  U 

the  downstream  A2  ^  (free  trailing  vortex  sheet)  in  the  para¬ 
meters  of  the  two  sides  (indicated  separately  by  "+"  and  "-") 
there  is  the  following  two  joined  conditions  [6,7]: 

P-  “  />■*•  (WA_-/U)l 

(/!.)-  -  (A.)+  -  0  /  ' 

Noting  that  A+  and  are  different  directions,  if  we 

let  =  _+ /i  +  t  and  AJ\.  =  A. /i  _>  then  when  the 

trailing  vortex  strength  is  equal  to  IJA/  and  the 
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direction  is  perpendicular  with  AJ\  then  it  is  parallel 
with  SA  .  Because  of  this,  generally  speaking,  the  comconents 
of  A  _  and  J\  always  show  outstanding  discontinuities  on 
the  free  trailing  vortex  surface.  Theoretically,  the  blade's 
downstream  must  always  be  able  to  show  a  free  trailing  vortex 
surface  for  all  of  the  non-free  vortex  flow  type  impellers. 

At  this  time,  it  is  only  necessary  to  accord  with  the  general 
Kutta-  >KYKOBCK  m  conditions  (abbreviated  as  the  K-J  con¬ 
ditions)  [14]  and  equation  (6)  in  order  to  derive  the  position 
and  flow  field  of  the  free  trailing  vortex  surface.  After  this, 
we  can  then  determine  2  d) .  Actually,  in 

order  to  satisfy  the  K-J  conditions,  we  can  make  the  starting 
line  of  the  free  trailing  vortex  surface  on  a  high  designated 
trailing  edge  stationary  point  position  along  the  blade. 

3)  On  the  Ag  surface.  It  should  satisfy  the  shock  discon¬ 
tinuity  conditions  [15] . 


(p<0-  —  OO+ 

(p» i)+ 

Or)-  —  Or)+ 

R-  -  R+ 


(7  A) 


From  this,  we  can  deduce  the  well-known  Rankine-Hugoniot  con¬ 
ditions. 


*  + 1, 

fp*\ 

\p~j 

l-i 

p-  i±j 
K-  1 

•-I 

(7B) 

\ 


Here,  the  subscripts  and  "+"  separately  indicate  the  aero¬ 
dynamic  parameters  of  the  shock's  front  and  back  sides. 
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Ill-  The  Variational  Principles  for  Three-Dimensional  Tran¬ 
sonic  Flow 


We  can  set  up  the  following  variational  principles  for  the 
above  mentioned  three-dimensional  flow.  The  solution  of  a 
three-dimensional  transonic  flow  with  embedded  shocks  in  a 
rotating  impeller  is  equivalent  to  the  extreme  value  nec¬ 
essary  conditions  of  the  following  function  J(«f  ): 


a]  —  o 


;(0)  -  jjj  F  (£,  |* ,  f*.  r)  I*  •  *  •  if  * 

+  jj  (*.)„,  *  iA-  j  j  (*"  -  *  -  A*)OA.)'  . 

—  jj  (♦+— 0-XM,)*  •  dA 


(8) 


In  the  equation 


il  j  ii 

See  fig.  1  for  the  A'  and  Ai'  surfaces.  The  symbols  and 

f  U  fa  /  U 

II  ||  II 

separately  indicate  the  parameters  of  the  A^  and  a"  surfaces, 

(V)  indicates  the  solution  domain,  and  d  indicates  the  position 
vector  of  the  interface  element  dA. 

[Proof]  When  we  take  the  variations  of  equation  (8)  and 
allow  each  position  of  the  unknown  discontinuities  to  be  free 
variations,  we  should  therefore  use  the  variations  with 
variable  domain  equation  (1-4)  in  the  Appendix  and  should  use 
stronger  conditions  £+=  #_  for  shock  surface  A  .  Because 

variables  and  sy  are  completely  arbitrary,  from 

£  J=0  we  can  thus  obtain  the  following  complete  set  of 
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extreme  value  necessary  conditions  : 

The  Euler  equation:  V*  ("p^l  )  =0  (which  continues  equa¬ 
tion  (1)  or  (1 ' ) )  . 

Group  of  natural  boundary  conditions: 

1)  On  A1  surface 


44.  -  (?.)„ 


2)  On  the  AJ  and  A'1  surfaces 

^  ^  U  /  u 

(44.)'  “  (44.)",  9"  —  G'  ■+>  a4>. 

rr‘rr.  This  it  can  he  inferred  that  their  tangential  fi 
rates  mast  ce  equal,  and  furthermore  from  equation  (3)  we 
.-mow  a.vat  4,4  and  all  aerodynamic  are  all  equal. 

3)  On  the  A_  ,  surface 

2 ,  d 

mu  -  mA.r  -o,  ?-r 


4)  On  the.  2  surface 


Cm.)-  -  054.)+,  $+ 


From  this  we  can  deduce  that  all  of  the  aerodynamic  para¬ 
meters  are  continuous. 

5)  On  the  A  shock  surface 
s 


(44.)-  ■  (44.)+ 

a* «■•).-( 

ow 
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On  the  Ag  surface,  besides  these  two  natural  boundary  condi¬ 
tions,  there  are  also  two  related  equations.  The  first  is  the 
above  used  strong  condition  3>_=  #  and  thus 

( a*  \ 

V  dV  )  [  dT  J+  wherein  the  tangent  flow  velocity  is  equal  to 


C^t)-  “  (/!,)+ 


(9b) 


The  other  related  equation  is  the  generalized  enthalpy  non¬ 
changing  conditions 


R-  -  R+ 


(9c) 


This  is  also  the  non-changing  relatively  stagnant  enthalpy 
★  ★ 

'V.-^wV 

Because  the  entropy  increases  in  transonic  flow  form  a 
direct  ratio  with  shock  strength  of  the  third  degree,  it  is 
only  necessary  for  the  Mach  number  before  the  shock  to  not 
be  much  higher  than  one.  Then  the  entropy  increase  after  a 
shock  will  be  very  small  and  the  changes  of  pQ  andfl^  can  be 

disregarded.  Therefore,  the  natural  boundary  conditions  of 
equation  (9a)  can,  along  with  the  conditions  of  equations 
(9b)  and  (9c) ,  represent  the  shock  discontinuities  conditions 
of  equation  (7A)  and  the  Rankine-Hugoniot  conditions  of 
equation  (7B) . 

Up  until  now,  we  have  derived  from  the  variational 
principle  »/(♦)  — o  the  potential  function  equation  (1*)  for 
transonic  three-diaensional  flow  with  shock  waves  in  the 
impeller  as  well  as  all  of  the  boundary  conditions  (including 
the  discontinuous  conditions  of  the  unknown  discontinuities). 
Moreover,  all  the  discontinuities  were  obtained  automatically 
using  the  forms  of  the  natural  boundary  conditions. 
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IV.  Generalized  Variational  Principle  for  Three-Dimensional 
Transonic  Flow  and  Its  Derivative  Group 


1.  We  can  also  establish  the  following  generalized  varia¬ 
tional  principle:  the  solution  for  the  three-dimensional 
transonic  flow  with  embedded  shocks  in  a  rotating  impeller 
is  equal  to  the  extreme  value  necessary  conditions  (here  let 
,  Ar>  A<?  .  Ay  ,  'p  and  be  independent  variables  of  the 

JT  Z 

following  functional  Jn(  <£  .  A  A?  ,  A  'p  and  /O  ) 

Ll  IT  Z 

7o(<P,  -If.  /!,.  A„  p ,  p)  —  jj  j  P  ■  r  •  Jr  •  J<p  •  Jg  +  j  j  •  JA 

(l’)  U\) 

~  J]  ( <P"  —  <t>-  A<t>J(pA.yjA 

-  [[(♦*  -  *-)(j>An)-4A  (10) 

(I) 


In  the  equation 


l  dr  r  d <p  9*  2 

11  "  *  (f  [*  “  ^  (^)l  ”4 


2/1./1,)} 


n 


[Proof]  Taking  the  variations  ofequation  (10) ,  and  using 

the  variations  with  variable  domain  of  formula  (1-5)  in  the 

Appendix  and  the  strong  conditions  of  $  ,  we  can 

derive  a  complete  set  of  extreme  value  necessary  conditions 

from  £  J  =0  as  follows: 
o 

The  Euler  equation  group.  Besides  including  the  four 
equations  in  (1)  and  (2) ,  it  also  has 

=1  (energy  equation) 

P-2T*"  (uniform  entropy  relation) 


b|l 


(A’  -43 
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They  are  equivalent  to  formula  (3) . 

The  natural  boundary  condition  group.  This  is  completely 
identical  to  the  previous  variational  principle. 

2.  Subgeneralized  variational  principle  group 

Following  the  conversion  method  found  in  reference  [8] , 
we  can  derive  a  family  of  subgeneralized  variational  principles 
from  this  generalized  variational  principle.  We  will  not  re¬ 
peat  the  individual  examples  here  but  only  point  out  that  among 
them  there  is  also  included  the  above  mentioned  variational 
principle  of  & J {  $>  )=0. 

V.  Conclusion 

In  the  above  section  we  established  a  family  of  varia¬ 
tional  principles  for  transonic  steady  relative  flow  with  em¬ 
bedded  shocks  in  a  rotating  impeller.  Its  special  feature  was 
taking  full  advantage  of  natural  boundary  conditions  and  all 
of  the  boundary  conditions  (including  the  discontinuous  condi¬ 
tions  on  the  various  unknown  discontinuities)  were  converted 
into  natural  ones.  This  provided  a  theoretical  basis  for  a 
new  computational  method  and  after  being  combined  with  the 
discontinuous  finite  element  [5]  as  well  as  the  special  hand¬ 
ling  of  the  supersonic  domain  we  were  able  to  use  the  finite 
element  method  to  automatically  and  distinctly  compute  the 
various  unknown  discontinuities  in  the  flow  field. 

Naturally,  the  results  of  this  article  are  likewise  applic¬ 
able  for  static  three-dimensional  cascades  as  well  as  various 
internal  flows  (for  example,  various  curved  flow  channels  such 
as  nozzles,  diffusers,  impeller  intake  and  exhaust  pipes)  and 
various  external  flow  [7]  problems. 

The  results  of  this  paper  can  be  viewed  as  extensions  of  the 
variational  principles  found  in  reference  [8]  and  as  related  in 


13 


sections  II  and  IV  to  three-dimensional  conditions. 
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Appendix  Analogy  Between  the  Functional  Variations  With  Variable 
Domain  and  the  Reynold ' s  Transport  Theorem 

Below  we  will  point  out  the  analogy  and  relationship  between 

variations  with  variable  domain  and  the  Reynold ' s  transport 

theorem.  Furthermore,  we  will  use  it  to  provide  a  simple  and 

convenient  method  for  solving  the  problem  of  variations  with 
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variable  domain.  Because  the  Reynold's  transport  theorem 
represents  the  relationship  of  the  control  and  closed  systems, 
we  can  use  a  mechanical  method  to  derive  it.  Therefore,  this 
method  is  more  directly  perceivable  than  most  methods  113] 
which  use  mapping  (conversion)  and  it  also  gives  variations 
with  variable  domain  a  certain  physical  (analogous)  signifi¬ 
cance. 

First,  Reynold's  transport  theorem  expresses  moving  liquids 
which  use  a  flow  velocity  w'  (12]  as: 

+  §  *•<*•">  0-0 
(V)  KV)  CCD 


In  this  equation,  (CS)  indicates  the  control  surface,  dV  is  the 
infinitesimal  volume,  dA  is  the  infinitesimal  surface  area,  the 
normal  direction  n  outside  the  (CS)  is  taken  as  regular,  and 
F=F(r,  <P  ,z,s)  is  an  arbitrary  function. 

We  can  compare  the  variations  of  the  area  to  the  position 
changes  of  the  closed  system  (V) .  The  position  changes  depend 
on  the  time  parameter  t.  Therefore,  for  the  arbitrary  variable 
52  ,  we  can  separately  write  the  following  formulas  for  its 
variations  with  fixed  domain  S  x  and  variations  with  variable 
domain  SQ  : 


8Q 


Naturally,  the  (CS)  moving  velocity  w'  can  be  indicated  by  its 
displacement  DS' 
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s 


<«*•»  *) 


i,  •  9r  +  >*ri<P  +  •  9m 


and 


**  -  * :4A  -  (''  s  +  ^tS* +  *'•  £)" 

After  using  t  to  multiply  each  item  in  equation  (1-1) ,  taking 
the  limit  as  t-*0,  and  taking  into  account  the  above  equations, 
we  obtain 

9  f  f f  f  •  dV  «■  [J|  9P  •  dV  +  <|  f{8*  .  a )4A 

(V)  icVi  (CJ> 

"•"  +  #f(£"  +  -5r^  +  x 

(CD 

Letting  J=  JJJ  FdV,  then  its  fixed  domain  and  variable 
(V) 

domain  variations  can  be  separately  expressed  as 

9j-9^r.4v 

icri  <v» 

After  substituting  in  the  above  equation,  we  can  obtain 

91-Sl  +  &  P(9f  .  n)dA  (I-J) 

(CD 


-JH 


9^jdA I 
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In  this  equation 


“  -  fjl  {£■  “  +  £  "■} " 


and 


*.  - 


<h»  ’ 


*.  - 


d* 


Introducing  vector  G 


c 


df 

a*. 


+  — 
a 


ar  •  J 
(*-)" 


Moreover,  after  applying  the  Gauss  theorem,  we  can  write 
-X 

S  J  as 


Sj 


dV  +  §  (C  •  <t)£*  .  dA 

<C»» 


After  substituting  in  equation  (1-2) ,  we  obtain 

tJ “  fJJ  (-S’  ”  v  * c)**  * ^  +  f  (c .  *•  + 

<e,,,  <CJ> 

Here,  S#  is  the  fixed  $  variation  on  (CS) .  In  order  to  make 
it  practical  and  convenient,  the  variation  along  the 

variable  domain  peripheral  interface  should  be  used  to  replace 
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Fig.  2  The  Variations  of  Function  ^  Along  the  Outside  Normal 
Direction 


Key:  1.  Fixed  CS  position 

2.  Peripheral  interface  position  after 
the  variable  domain 

From  fig.  2  we  know  that  6$  =  £  +^  =  Sx  +  £5^, 

and  after  substituting  in  formula  (1-3) ,  we  obtain  a  general 
formula  for  variations  with  variable  domain 

«-  JJJ  {JZL  .cy».Jv  +  £  {e. .  «•  +  (i-*) 

'  (CS) 

When  F  contains  many  unknown  function  =1,2,3...),  then 


"  -  HI  2  [(£■  -  » •  o’)"*]  *  f  2  {<*  •  *•' 


(1-5) 
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VARIATIONAL  PRINCIPLES  AND  GENERALIZED  VARIATIONAL 
PRINCIPLE  FOR  THREE-DIMENSIONAL  TRANSONIC  FLOW 
WITH  SHOCK  WAVES  IN  A  ROTATING  TURBO-IMPELLER 


Liu  Gaolian 

(Shanghai  Inttitute  of  Meek.  Sngintering) 


Abstract 

In  this  article  a  family  of  variational  principles  for  3D  transonic  steady  relative 
flow  with  embedded  shocks  in  a  turbomachine  impeller  of  axial,  radial  or  mixed  flow 
type  is  developed.  Its  special  feature  is  to  take  full  advantage  of  natural  boundary  con¬ 
ditions  and  “artificial  interfaces,”  so  as  to  facilitate  the  handling  of  various  complex 
boundary  conditions.  As  a  result,  all  boundary  conditions  in  the  problem  under 
study  have  been  converted  into  natural  ones  Moreover,  it  is  also  shown  that  by  taking 
variations  of  the  position  of  unknown  flow  discontinuities  (such  as  shocks,  free 
trailing  vortex  sheets  all  matching  conditions  across  these  discontinuities,  including  the 
well-known  Rankine-Hugoniot  shock  relation,  can  be  derived  from  the  variational  prin¬ 
ciples  as  natural  interface-conditions. 

This  article  is  primarily  intended  to  provide,  in  conjunction  with  the  discon¬ 
tinuous  finite  element1’1,  a  theoretical  basis  for  developing  a  novel  computational 
method,  which  allows  all  flow  discontinuities  to  be  singled  out  automatically  and 
distinctly  (without  smearing).  Owing  to  the  assumption  of  the  existence  of  a 
velocity  potential,  the  applicability  of  the  present  theory  is  limited  to  cases  where 
Mach  numbers  before  the  shock  do  not  differ  greatly  from  one.  The  variational  princi¬ 
ples  presented  herein  are  some  extension  of  those  developed  previously  in  (8J. 

In  the  appendix,  an  analogy  between  the  functional  variations  with  variable 
domain  and  the  Reynolds’  transport  theorem  is  pointed  out  and  used  as  a  basis  for 
deriving  a  general  formula  for  such  variations. 
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CALCULATION  FOR  THREE-DIMENSIONAL  TURBULENT  BOUNDARY  LAYER  ON 
A  YAWED  WING  IN  COMPRESSIBLE  FLOW 

by  Zhang  Guofu  „  ‘ 

(Nanjing  Institute  of  Aeronautics) 

Abstract 

This  paper  starts  from  the  equations  of  motion  in  gener¬ 
alized  curvilinear  coordinates  and  the  model  of  Head's  two- 
dimensional  entrainment  theory,  and  derives  a  set  of  equa¬ 
tions  along  surface  curve  on  normal  profile  of  a  yawed  wing 
for  a  compressible  three-dimensional  turbulent  boundary  layer. 
After  a  typical  form  of  a  set  of  ordinary  differential  equa¬ 
tions  with  initial  values  is  obtained  by  matrix  transform¬ 
ation  we  solved  these  equations  by  numerical  computation. 

Thompson  and  MacDonald's  computational  results  on  flow 
around  an  attachment  line  in  the  leading  edge  of  a  yawed  wing 
have  been  used  as  the  initial  values  to  solve  the  set  of 
ordinary  differential  equations. 

The  method  in  this  paper  can  be  used  to  calculate  the 
separation  position  of  boundary  layer,  various  boundary  layer 
thicknesses  and  wall  stress,  and  describe  the  difference 
between  streamlines  outside  the  boundary  layer  and  streamlines 
on  the  wall. 

Typical  computations  show  that  the  present  method  agrees 
well  with  experimental  results. 
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Primary  Symbols 


u 

c 

local  composite  velocity  on  boundary  layer 

u,v,w 

velocity  in  £  ,n,§  directions 

VV1 

component  velocity  of  Uc  in  x,y  directions 

CM 

*H 

lateral  stream  angle  (fig.  1) 
local  stream  line  angle  (fig.  2) 

=  tge 

shearing  stress  along  *b  and  q  directions 

^  w 

wall  shearing  stress 

*Vl  ,"5at2 

components  in  the  \  and  n  directions 

w  1 

Cf 

T'  Ah.fi  U2)  ,  friction  coeffient 
w/  e  e 

cfl 

=cfcos 8,  the  friction  coefficient  in  the  main  stream 

direction  outside  the  boundary  layer 

Cf  2 

=  cfcos£,  the  friction  coefficient  in  the  lateral 
stream  direction 

V 

A 

air  motion  viscosity  coefficient 
wing's  angle  of  yaw 

c 

normal  wing  profile  chord  length 

thickness  of  boundary  layer 

=  \S  (l-jjj^y  )  displacement  thicknesses  in 
Jo  e  e 

and  f)  directions 

d  4, 

s2 

cS  — 

=  -  \  4»v  displacement  thicknesses  in  d 

Jo  U 

and  directions 

=  ^ (1-  ^  momentum  loss  thicknesses 

Jo  u«  *cse 

:  in 

d£  and  ^  directions 

CM 

<D^ 

=  -  ^  ^VV-  momentum  loss  thicknesses  in 

Jo  P.u. 

directions 

d  and 
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★ 

c 


H 


H 


o 


mixed  momentum 
loss  thicknesses 


-  (Vu)/v. 


=  ’  (3/3/dx) 

Vo£  sin^/1 


attachment  line  in  air  stream 
with  similar  parameters 


=  s 
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compressibility  analogy  of 
d  and  H 


h  =  (  S  -  8  ,>/  e 

1  1  ii 


Footnotes 

a.l  attachment  line  value 

e  potential  flow  value  of  boundary  layer 

w  wall  value 

A  lateral  line  on  top  of  a  symbol  indicates  a  time  average 
value  and  one  with  a  left-falling  stroke  indicates  a 
pulsating  value. 

I .  Foreward 

In  recent  years,  due  to  developments  in  aviation, 
especially  those  in  transonic  flight  technology,  aero¬ 
dynamics  is  daily  confronted  with  the  serious  challenge 
of  separated  flow.  For  example,  in  mechanism,  the  solution 
of  the  problem  of  transonic  buffeting  being  totally  differ¬ 
ent  from  low  speed  buffeting  always  possesses  a  decisive 
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effect  on  the  economic  performance  of  the  transporter 
(buffeting  boundary)  or  on  the  mobility  performance  of  a 
fighter  (buffeting  performance) . 

When  the  boundary  layer  theory  is  used  to  handle  the  pro¬ 
blem  of  three-dimensional  wing  separation,  even  the  simplest 
calculation  of  an  unlimited  yawed  wing  boundary  layer  is  much 
more  difficult  than  for  a  profile.  There  are  two  main  reasons 
for  this:  (1)  When  considering  the  viscosity  of  the  flow 
around  an  unrestrained  yawed  wing,  the  calculation  results  for 
a  normal  profile  can  usually  be  directly  converted  to  a  yawed 
wing  in  accordance  with  the  law  of  cosines.  Yet,  when  handling 
the  related  viscosity  problems,  in  principle  the  law  of 
cosines  is  not  established.  This  is  because  the  S  shaped 
curved  stream  line  outside  the  boundary  layer  not  only  causes 
the  inside  of  the  boundary  layer  to  have  the  same  velocity  as 
the  main  stream  direction  outside  the  boundary  layer  but  the 
inside  also  has  the  same  vertical  velocity  as  the  main  stream 
(fig.  1) . 

(4 ) 


Fig.  1 

Key:  1.  Wall  stream  line 

2.  Boundary  layer  stream  line 

3.  Wall  stream  line 

4.  Boundary  layer  stream  line 

5.  Separating  line 
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Therefore,  the  stream  within  the  boundary  layer  is  three- 
dimensional.  Concurrently,  the  wall  stream  line  is  not  identical 
to  the  boundary  layer  stream  line  and  the  included  angle  between 
them  is /S-  Any  one  point  on  the  wall  takes  the  coordinates  i|  , H , 

3  4?  is  the  stream  line  direction  along  the  boundary  layer, 

H  is  vertical  to  along  the  wall  and  is  upward  vertical 

to  the  wall.  Fig.  1  draws  the  main  stream  direction  velocity 

u(  £  )  and  the  lateral  stream  direction  velocity  v(£)  with  the 

corresponding  wall  friction  stresses  of  and  •  The 

wl  w2 

total  friction  stress  is  thus  =  '£',+  1?  _.  The  wall  stream 

w  wl  w2 

line  direction  is  defined  as  identical  to  the  'f  direction. 

w 

Generally,  Redeker  [1]  points  out,  for  a  wing  with  a  relatively 

large  angle  of  sweepback  (  /i  ^  30°) ,  there  is  separation  when 

the  wall  stream  line  is  parallel  to  the  trailing  edge.  Although, 

a-t  this  time,  the  wall  stress  is  not  equal  to  2ero,  yet  for  a 

relatively  small  wing,  the  occurence  of  separation  is  often 

already  7Tw=0  before  the  wall  stream  line  is  parallel  to  the 

trailing  edge.  This  is  analogous  to  a  two-dimensional  situation. 

(2)  When  the  air  stream  flows  along  the  leading  edge  of  the 

yawed  wing,  two  branches  are  formed  on  the  attachment  line,  one 

branch  goes  up  while  the  other  goes  down.  This  phenomenon  is 

analogous  to  the  air  stream  flowing  gz? lually  passed  and  closer 

to  a  drawing  in  flat  boundary  layer.  Based  on  the  test  results 

of  Cumpsty  and  Head  [2] ,  when  ^  100  is  on  the  attachment 
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line,  there  is  a  turbulent  boundary  layer.  When  modern  large 
scale  transporters  are  cruising,  they  generally  exceed  this 
value  and  therefore  we  can  consider  that  there  is  a  turbulent 
boundary  layer  on  the  attachment  line.  By  utilizing  a  series 
of  dominant  equations  for  the  turbulent  boundary  layer  and 
integrating  the  special  characteristics  of  the  stream  on  the 
attachment  line,  we  can  calculate  the  air  stream  parameters  on 
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the  attachment  line.  Crumpty  and  Head  [*2]  initially  calculated 
uncompressible  results  but  compressible  and  more  practical 
results  were  subsequently  presented  by  Thompson  [3] . 

Due  to  the  difficulties  of  test  research  and  analysis  of 
the  above  mentioned  theory,  progress  in  this. area  has  been 
slow.  Not  until  the  1970‘s  did  we  begin  seeing  research  papers 
on  large  swept  back  wings  with  compressible  turbulent  boundary 
layers.  Yet,  some  papers  only  present  formula  results  and 
others  are  merely  limited  to  discussing  specialized  problems. 
They  lack  systematic  formula  demonstration  and  calculation 
methods.  In  addition,  many  errors  have  been  found  in  the  formulas 
published.  In  view  of  the  above  mentioned  facts,  this  article 
starts  from  the  equations  of  motion  in  generalized  curvilinear 
coordinates  and  the  model  of  Head's  two-dimensional  entrainment 
theory  and  derives  a  set  of  dominant  equations  for  the  stream 
line  in  the  curvilinear  coordinates  along  the  boundary  layer. 

The  results  are  basically  the  same  as  those  found  in  reference 
[1] .  Redeker  [1]  used  repeated  iteration  to  solve  the  above 
mentioned  set  of  equations.  Afterwards,  Thompson's  [3]  research 
showed  that  the  iterative  method  could  have  divergent  problems 
when  applied  to  the  leading  edge  of  a  relatively  large  wing 
with  an  accelerating  pressure  gradient.  Moreover,  convergence 
was  very  slow  near  the  separation  area.  We  think  that, 
practically  speaking,  to  carry  out  numerical  calculations  of 
the  stream  line  along  the  outside  of  the  boundary  layer  is 
not  as  convenient  as  along  the  normal  profile  curve  and  there¬ 
fore,  in  this  paper,  we  have  carried  out  variable  substitution 
and  matrix  operation  so  as  to  obtain  a  typical  form  of  a  set  of 
ordinary  differential  equations  with  initial  values  along  the 
normal  profile  curve.  The  velocity  and  wall  friction  stress 
formula  selected  in  this  article;  used  the  conditions  on  the 
attachment  line  given  in  reference  [3]  as  initial  values; 
used  numerical  computations  for  their  solution. 
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As  regards  a  correctly  designed  transonic  wing,  wing  root 
and  wing  tip  modifications  should  be  carried  out  so  as  to 
maintain  as  best  as  possible  the  uniformity  of  their  isobaric 
lines  and  equal  percentage  chord  lines.  This  will  cause  the 
greater  area  on  the  transonic  swept  back  wing  to  possess 
unrestrained  yawed  wing  qualities.  Therefore,  the  boundary 
layer  results  from  research  on  the  unrestrained  yawed  wing 
will  have  real  significance. 

II.  Set  of  Dominant  Equations  For  a  Turbulent  Boundary  Layer 
With  Curvilinear  Coordinates 


We  can  rewrite  the  equation  of  single  quality  gas  motion 
in  generalized  curvilinear  coordinates  presented  in  reference 
[4]  for  situations  with  applied  force  (including  pressure  and 
tangential  force)  as 


_L[iL(8T\  §T] 

h,  i it  Va iJ  0<7,J 


- —  +.  p 

Ph<  8* 


(1) 


In  the  equation,  qi  is  the  generalized  curvilinear  coordinate, 
q^  is  the  derivative  of  q^  for  time  t,  F.^  is  the  tangential 
force  working  in  the  i  direction  component  on  the  single 
quality  gas,  T  represents  the  kinetic  energy  of  the  single 
quality  gas,  and  h^  is  the  Lame  coefficient  along  the  coordin¬ 
ate  q^  Then 


T 


^  j 

~~  2  ?i*  1**  4u  4*) 

2  i -i 


(2) 


Taking  the  orthogonal  curvilinear  coordinates  £  ,  <7  and  £ 
shown  in  fig.  1  and  letting  hj=l,  from  formula  (2)  we  then 

obtain 
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(3) 


9T  _  ul  dki  t  »*  6J,, 
Qq,  h\  dqi  A,  9f, 


i.  ( H\  _  y>  “i.  Ph  +  h  j9“A 

dt\diJ  ^ih,\'dq,  '  dqj 


hi 


du, 


0) 


In  formula  (1),  can  be  shown  by  the  deformation  angle 
velocity  of  the  curvilinear  coordinates  [4]  and  using  Prandtl 1 s 
boundary  layer  assumption  we  can  obtain 


Pressure  strength  p.  can  also  be  expressed  by  the  pressure 
strength  arithmetic  mean  p  [5]  along  the  three  coordinate 
directions . 

In  turbulent  motion,  we  can  regard  the  real  gas  flow  para¬ 
meters  as  being  formed  from  the  superposition  of  the  time 
average  and  pulse  values,  and  based  on  the  time  averaging 
principle: 

z  -  2  +  r  . 

_  _  \  (6) 

ST  —  (5  +  S’X  f  +  T*)  -  ST  +  YF  f 

If  we  substitute  formulas  (3) ,  (4)  and  (5)  into  formula  (1) , 
carry  out  time  averaging  for  formula  (6) ,  and  then  carry  out 
simplification  based  on  Prandtl ' s  boundary  layer  assumption, 
we  will  finally  obtain  a  set  of  equations  of  time  average 
motion  in  curvilinear  coordinates  for  a  turbulent  boundary  layer 
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%  direction  momentum  equation 


,  —  du  .  ,  —  da  ,  ,  ,  —  du  ,  — - 
—  +  htf>v  — —  4-  —  +•  XT 

Of  3*1  O?  Oil 


+  AiAj 


(7) 


0  direction  momentum  equation 


(8) 


In  the  formulas 


ri  —  /»(05/9f)  —  (p<f)V  (9) 

tx  “  p(0S/ Of)  —  «V  (10) 


The  first  item  on  the  right  side  of  the  above  two  formulas 
is  the  tangential  stress  caused  by  viscosity,  and  the  second 
item  is  the  tangential  stress  (called  Reynold's  stress)  pro¬ 
duced  because  of  the  velocity  pulse  when  there  is  turbulence. 
Physically,  the  above  results  are  identical  to  the  conclusion 
when  there  are  rectangular  coordinates. 


From  the  above  similar  methods  we  carry  out  time  averaging 
for  the  continuous  equations  [4]  of  the  generalized  curvilin¬ 
ear  coordinates  and  obtain 

+  +  (ll) 

If  we  calculate  the  %  coordinate  based  on  the  velocity 
potential,  the  d  |  =d9<  .  If  along  the  %  coordinate 

dsssh^dfc  ,  we  obtain 
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At  -  VJU, 


(12) 


;j 

Then,  along  the  boundary  layer,  formulas  (7)  and  (8)  can  be 
simplified  as 


_ _L  it] 

A,  d§  hx  d§ 

u.  du.  _  i  a? 

***  Aj  di)  A]  di)  J 


03) 


Fig.  2 

Key:  1.  Wall  streamline 

2.  Boundary  layer  streamline 

As  shown  in  fig.  2,  the  normal  profile  curve  along  the  wing 
takes  coordinate  x,  x  begins  from  the  attachment  line,  and 
takes  coordinate  y  along  the  span.  The  included  angle  between 
the  wall  stream  line  and  boundary  stream  line  is  &  (lateral 


stream  angle) ,  and  the  included  angle  between  the  4  coordinate 
and  x  coordinate  is  £  (local  stream  line  ancle) .  The  component 
velocity  of  in  the  x  and  y  directions  is  u^  and  v^.  Now  we 
will  consider  the  following  relations: 


(1)  Because  it  is  an  unrestrained  yawed  wing  and  the  air 
stream  parameter  along  the  y  direction  is  unchanging,  we  can 
obtain 


a 


-tga 


±_d 
hx  a? 


-K 


_a 

a. 


.a  A . 

«i  ds 


(H) 


(2)  Based  on  the  formula  defined  by 


and 


we  obtain 


0,2+-a,-02I  (15) 


(3)  when  we  consider  compressibility,  there  is  a  related 
formula  of  the  stream  line  along  the  boundary  layer 


JL  _ >f2  1  dU. 

P.  ds  ’V,  ds 


(  16) 


If  we  take  the  thickness  integral  of  each  item  in  formulas 
(7) ,  (8)  and  (11)  along  the  boundary  layer  and  substitute  them 
into  formulas  (12-16) ,  after  operation  we  obtain: 


The  %  direction  momentum  equation 


WiL+dii  +  M»)  ^7T~Kdxs 

ds  v.  ds  v  >  ds  v,  ds 


+  -i-  fi  (0U  -  eu)  -  -£■ 

A,  ds  P,Vl 


(17) 
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The  1  direction  momentum  equation 


dd:1 

ds 


PJL  (2  -  Ml)  +  0U  -1  dltJ  -K  &L 

U.  ds  A,  ds  ds 


K  dU± 
Ut  ds 


■  Mi  ea- (ea  +  «„  +  *,) 
U  0  Os 


(18) 


We  took  Head ' sAtwo-dimensional  entrainment  theory  and  ex¬ 
tended  it  to  a  compressible  three-dimensional  orthogonal 
curvilinear  coordinate  system.  After  the  unit  area  on  the 
boundary  surface  was  put  into  the  gas  stream  of  the  boundary 
layer,  we  used  the  entrainment  function  FtH^)  to  express  it 


1  /  do_  +  dt_ 
u.\h t  9{  A,  dr, 


w.)  -  F(H.) 


(19) 


We  took  the  thickness  integrals  of  each  item  of  the  contin¬ 
uous  equation  (11)  along  the  boundary  layer  and  applied  the 
dif f erentio-integral  theorem  to  obtain 


_L  JL  dC  j_  _a_  [‘heitjc 
A, A,  df  J»  U.  hxhx  a,  J •  u. 


£t  9JL 

a,  a? 


hi 


dS  \ 


(20) 


If  the  right  side  of  the  above  formula  is  substituted  with 
formula  (19) ,  the  left  side  undergoes  differential  operation 
and  the  independent  variable  position  is  converted  to  s,  we 
finally  obtain  the  entrainment  equation 


-  Mix*  -•»)  +  -T  ¥* («-«»)  +  *  90 


ds  u.  Ot  '  "N  A,  ds 

f  (ho 

V,  ds 


ds 


(21) 


In  the  formula,  the  entrainment  function  F(H^)  uses  the  two 
dimensional  empirical  formula 


F(H,)  -  0.0299(H,  -  3)^  "7 
Hi  -  1.535(H  -  0.7)“*-"»  4-  3.3 


(22) 

(23) 


In  order  for  the  various  boundary  layer  thicknesses  in  the 
above  formulas  to  be  able  to  have  intggration,  we  had  to  know 
the  stream  velocity  and  lateral  velocity  types  of  the  turbul¬ 
ent  boundary  layers.  We  used  the  velocity  type  formula  recom¬ 
mended  in  reference  [3]  which  calculated  the  compressibility 
effect : 


The  lateral  velocity  type  based  on  multiplication  rate  form 


U. 


(24) 


The  lateral  velocity  type  based  on  an  extension  of  the  non-com- 
pressible  stream  formula 


V. 


(i  -JL)± 
\  zju. 


>gA 


(25) 


In  the  formula 


(fl~  l)/2 
ff  -  7,(H,)«  1 

H-y, (3,  M.)«J 


(26) 

(27) 


From  formulas  (24)  and  (25) ,  Green  integrated  the  various 
boundary  layer  thicknesses  as^J,  and  n  functions.  After 

trying  formula  (26)  we  obtain 


9m  -  -  HMAWye u  i 

«.-«§ *f,(ff)*n;  i  («) 
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In  the  formula 


1 


U  “  — 2/[(H  —  l)(ff  +  2)] 
h  -  (14/?  +•  30)/[(H  +  2)(H  +■  3)(/7  +  5)] 

/3-  -lfiff/[(ff-  1)(/?  +  3)(/?  4-  5)J  (29) 

U - 24/[(/?  -  1)(/?  +  2)(/?  +  3)(H  +  4)] 

ti-U-U  ' 

In  order  for  application  to  be  convenient  the  stream  line 
coordinate  s  on  the  outer  part  of  the  boundary  layer  was  trans¬ 
formed  into  the  curvilinear  coordinate  x  along  the  normal 
profile.  We  used  the  divergence  formula  of  the  main  stream 
velocity  on  the  boundary  layer  in  the  (4  >  H  ,  £  )  and 
(x,y,z)  coordinates. 


div  (i/,) - - — £LE£ 


;4|)  __  dut  di>i 

6*  dy 


and  obtained 


J_  .  JL  (§ax  +  _ L_  ( 

*1*.  df  U.\Ox  dy  )  v*x  d$  1 

For  the  unrestrained  yawed  wing  (  9  9y)=0,  and 


d  _  ix  d  M  d 
df  4,  dx  U,  dx 

By  utilizing  the  above  mentioned  relation  we  then  simplified 
formula  (30)  into 

_L  -  l  3*,  _  rMn  _  l  W  dU. 

Ma  e?  fi,  d,  U\  dx  Ul «,  dx  v  ' 


We  substituted  formulas  (28),  (31)  and  (32)  into  formulas 
(17) ,  (18)  and  (21)  and  finally  obtained  a  set  of  dominant 
equations 
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Ku  K\i  Ku\ 

K- «» *•)  ir 

Ks,  Kj  d(g  _  g  ) 


In  the  formula 


^  _  £1  f  t g  p  +  £i  tg^  H 

V,  U,%*  V.  *  dR  JHt  1 

Ku - jftA, 

U  0 

K»  -  ?/■"'  -  g  '•*'-*  (g  «S 

-a.g^  '(i) 

t/,  <//// 

Ku-^fAi  -  -LjfUdu'gt 
v  t  U  g 

K  -  tg/?  -  -fl  tg> £ 

V.  dH  uHx  U.'dH  dHt 

-Z"'* 

a  fAi 

V  9 

K«  -  a  +  a  t g(y  -  (l  -iE 

£/,  t/,  £«  £W, 

F‘  "  ^  ~  §*  [(»  +  2  -  Mi)  +  ^  (1  _  * 

+  KUW  a]  -fflt 

“  tg/?¥  +  $*“['*'  •  /.(mi  - 2  S) 

+  ff(W+l  +  /«tgJ£(l  -Mi)  )] 

*  -  F(Ht)  +  a  tfu  [h,  (Mi  ~^)  +  MitgV/ij 
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In  the  formulas,  dH/dH^,  H  and  H  can  be  expressed  by 
formula  (27)  as  the  function  of  .  c^  can  be  expressed  by 
the  empirical  formula  recommended  in  reference  [3]. 

The  set  of  equations  in  (33),  through  matrix  transformation, 
obtained 

M„  Mu  Mu\  /FA 

Mn  Ma  Ma  IF,)  (34) 

Mn  Ms  M \Fj/ 


In  the  formula 

Mn  “  K-aKtiy  Mu  “  XuA.'s 

Mu  “  KuKa  KuKu,  Mjj  *  ^a^ji  — 

Ma  -  k,iKb  -  a:,^,  m,  - 

Ms  "*  ^ ,  Ms  ■*  ^ijKh  KuKa 

MB  ■*  KnKji  —  KuKu 

1*1  -  Kh/CsA:*  -  KuKaKa  -  KuKnKm  +  K„KaK„ 
t  KoKnK„  —  KaKaKa 

Formula  (34)  is  a  set  of  closed  ordinary  differential 
equations  with  initial  values.  These  equations  can  be  solved 
by  the  Gill  or  R-K  methods. 

III.  The  Flow  on  the  Attachment  Line 

The  numerical  method  used  by  Thompson  [3]  to  solve  the  set 

of  dominant  equations  for  the  turbulent  boundary  layer  on  the 

attachment  line  expressed  R  g  ,  n  and  Z.  as  the  curvilinear 

*  11  1 

or  analytic  formula  of  C  and  M  . .  The  formula  used  was 

a  •  i 
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“  (1  +  2n)ln 


05) 


and  this  could  find  0^,  tgB  and  -  £  on  the  attachment 

line.  Because  u^=tgB=0,  that  is  )k)  =0  on  the  attachment  line, 
therefore  with  the  set  of  equations  in  (34)  there  is  singular¬ 
ity  on  the  attachment  line.  Because  of  this,  calculations 
should  begin  from  the  downstream  ^  x  area  on  the  attachment 
line.  In  light  of  the  fact  that  the  flow  on  the  attachment 
line  is  symmetrical,  then 

tgj?  -  (ddjdx)  -  d(a  -  ax)/dx  -  o 

Taking  B  tg^  and  S  -  8  ^  in  the  neighborhood  of  the 

attachment  line  in  accordance  with  Tylor's  expansion  in 
series,  we  retained  the  first  non-zero  item  and  obtained  the 
initial  value  in  the  4  x  area 

(#u)a«  “•  (tg^)a»  Ax  •  ; 

\  lx  /..i 

(*  “  *)*.  -  («  -  «.)*, 

The  calculations  in  this  paper  which  took  ^x=0.005cn, 
attained  satisfactory  results. 

IV.  The  Process  of  Numerical  Calculations  and  Calculation 
Examples 

1.  Based  on  the  potential  flow  theory  we  calculated  velocity 
distribution  u^  and  span  velocity  v^  outside  the  boundary  layer 

on  the  normal  profile  of  a  yawed  wing. 

2.  From  the  known  M  .  and  C*  values  along  the  attachment 

0*1 
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line  we  calculated  (  0  )  . ,  (dtg/j/d x  )  .  and  (  -  S  i  i 

i  a  3  •  i-  3..JL  ^  l  a .  1 

and  then  from  formula  (36)  we  calculated  the  above  mentioned 
values  of  the  downstream  A  x  area  on  the  attachment  line  so 
as  to  solve  the  initial  values  of  the  set  of  equations. 


3.  We  used  the  Gill  (or  R-K)  method  to  carry  out  a  numerical 

solution  for  the  set  of  equations  in  (34)  and  finally  we  found 
the  £  -  £“  as  well  as  £  and  the  wall  friction 

coefficient  values  on  each  nodal  point  along  the  X  coordinate. 

4.  When  we  calculated  £ =  90°  (fig.  2)  there  was  separa¬ 
tion  . 


Fig.  3  is  Cumpsty  and  Head's  experimental  measurement 
results  for  the  boundary  layer  of  a  yawed  wing  with  a  sweep- 
back  angle  of  60°. 

(*)  <h)  (e> 


(5)  *.0  WtoT-ii  - ftedefer  itJUg*, 

(6)  *'  "  -I.*-  io*. 

Fig.  3 

Key:  1.  Separation 

2.  Separation 

3.  Separation 

4.  Separation 
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Fig.  3  (continued) 

5.  #,  0  are  the  measurement  results  of  two  different 

types  of  boundary  layers, -  indicates  the  cal¬ 
culation  results  of  this  paper,  -  indicates 

Redeker's  calculation  results 

6.  Re=  —  =1.4*10  ,  trailing  edge  coordinate 
trailing  edge  coordinate  x=1.2c 

n 

Fig.  (3a)  is  a  sketch  of  an  experimental  installation,  fig.  (3b) 
is  a  normal  wing  profile,  fig.  (3c)  is  the  measured  velocity 
distribution  on  the  wing  surface,  figs.  (3d),  (3e)  and  ( 3 f )  com¬ 

pare  the  values  and  measurement  results  of  Q  ^ ,  H  and  /3  along 

x  presented  in  this  paper  with  the  calculation  results  given  by 
Redeker  [1],  fig.  (3g)  is  this  paper's  calculated  distribution 
value  along  x  of  the  wall  stream  friction  coefficient.  We  can 
see  from  this  figure  that  the  calculation  results  of  this 
paper  agree  well  with  the  experiment  results.  When  compared 
with  Redeker's  calculation  results,  the  results  given  in  this 
paper  for  the  area  near  the  front  of  the  separation  line  are 
even  more  identical  to  the  experiment  results.  This  paper  cal¬ 
culated  the  separating  point  coordinate  X=1.004cn  (a  little 

forward  as  compared  to  Redeker's  results)  ,  ,#=33 . 63° ,  h=1.755, 

0 ,  =c  *  8.412*10  .  In  the  calculation  examples,  a  laver 

11  n 

stream  recovery  area  did  not  appear  (reference  [10]  mentions 
the  possibility  of  a  recovery  area  appearing  under  a  very  strong 
leading  edge  accelerating  pressure  gradient) .  This  program  was 
computed  on  a  TQ-16  in  approximately  25  minutes  (including 
101  dot  output) . 

References 

[1]  Redeker,  G.Z.,  Z.  Flugwiss,  10(1973). 

[2]  Cumpsty,  N.A.,  Head,  M.R.,  Part  I, II,  Aeron.  Quart. 

18 (1967)  . 


37 


[31 

■4] 


Thompson,  MacDonald,  A.G.J.,  ARC  CP,  1307(1974). 

Ke  Qin,  H.E.,  Ji  Bieli,  N.A. ,  Luo  Si,  Theoretical 
Hydromechanics.  Vol .  1,  Book  1,  32-45. 


[5]  Sjh  LUh-h,  Tvp6y.ieHTH0e  reqeiuie  *HjKOCTeH  h  raioa,  120—126. 

[6]  Lu  Hanjiang,  The  Attached  Layer  Theory,  152-155. 

[7]  Heal,  M.R.,  ARC  R.&  M. ,  3152(1958). 

[8]  Green,  J.E.,  J.F.M.,  31,  Part  IV(1968) ,  753-778. 

[9]  Cumpsty,  N.A. ,  Head,  M.R.,  Part  IV,  Aeron.  Quart,  21(1970). 

[10]  Thompson,  B.C-.J.,  ARC  CP  1308  (1974). 


CALCULATION  FOR  THREE-D1MENTIONAL  TURBULENT 
BOUNDARY  LAYER  ON  A  YAWED  WING  IN 
COMPRESSIBLE  FLOW 

Zhang  Gaofa 

Abstract 

This  paper  starts  from  the  equations  of  motion  in  generalised  enrvilinear 
coordinates  and  the  model  of  Head’s  two-dimensional  entrainment  theory,  and  derives 
a  set  of  equations  along  surface  curve  on  normal  profile  of  a  yawed  wing  for  a 
compressible  three-dimensional  turbulent  boundary  layer.  If  there  is  a  strong 
favourable  pressure  gradient  in  the  neighbourhood  of  leading  edge  of  the  wing,  diver¬ 
gence  will  arise  in  eonvertional  methods.  Therefore  a  typical  form  of  a  set  of 
ordinary  differential  equations  with  initial  values  is  obtained  by  matrix  transformation 
in  this  paper.  These  equations  have  been  solved  by  numerical  computation. 

Thompson  and  MacDonald’s  computational  results  on  flow  around  an  attachment 
line  in  the  leading  edge  of  a  yawed  wing  have  been  used  as  the  initial  values  to  solve 
the  set  of  ordinary  differential  equations. 

The  present  method  can  be  used  to  calculate  the  separation  position  of  boundary 
layer,  various  boundary  layer  thicknesses  and  wall  stress,  and  describe  the  difference 
between  streamlines  outside  boundary  layer  and  streamlines  on  the  wall. 

Typical  computational  shows  that  the  present  method  egress  well  with 
experiment  results. 
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ON  THE  RAYLEIGH- TAYLOR  INSTABILITY  OF  A  CURRENT  SHEET  IN  AN 
ELECTROMAGNETICALLY  DRIVEN  SHOCK  WAVE  APPARATUS 

by  Xu  Fu 

(Institute  of  Mechanics,  Academia  Sinica) 

1.  Current  in  Shock  Wave  Apparatus 

At  present,  we  will  consider  the  electromagnetically 
driven  shock  wave  apparatus  as  shown  in  fig.  1.  After  dis¬ 
charge,  a  very  strong  shock  wave  is  propagated  towards  the 
right,  the  rear  surface  is  a  current  layer  and  it  separates 
the  plasma  and  magnetic  field.  In  1959,  Wright  and  Black  [1] 
researched  in  detail  the  current  in  an  electromagnetically 
driven  shock  wave  apparatus  during  the  initial  period  of 
discharge.  After  assuming  that  current  i  and  time  t  form  a  dir¬ 
ect  ratio,  they  found  a  similar  solution.  At  the  same  time,  the 
current  layer  used  iso-acceleration  to  move  towards  the  right 
side  and  it  could  be  predicted  that  this  would  produce  Ravleigh- 
Taylor  instability.  Because  acceleration  a  was  large,  the  in¬ 
stability  growth  rate  w  was  also  very  large. 


Fig.  1 
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Fig.  1 


Key:  1.  Electrode 

2.  Magnetic  field 

3.  Electrode 

4.  Current  layer 

5.  Plasma 

6 .  Shock  wave 


We  let  the  current  sheet  use  iso-acceleration  a  to  move 
towards  the  right  (see  fig.  2) 


Fig.  2 

Coordinate  xQ(t)  and  velocity  vQ(t)  of  the  curr  ent  sheet  are 
separately 


*o(0  4-  bt 

fa (/)  —  at  +  b 


0) 


Given  that  the  plasma  is  compressible,  non-viscous  and  non¬ 
heat  conducting,  the  electric  conduction  rate  is  (5  -  CO  , 
and  the  magnetic  conduction  rate  is  14=1.  Using  Gauss  units, 
then  the  fundamental  equations  are: 


—  —  v  x  (V  x  //) 
at 

v  ■  H—  0 

- vp-^Hx  ( vxff) 

dp 

-ft—  +  ( V  •  V)p  4-  pV  •  V  —  0 

+  ( V  •  V),  -  ,  -i  {  &-  +  ( V  .  v),} 


(2) 
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The  kinematic  bound  conditions  on  the  current  sheet  are: 

On  x=x  (t) 
o 

V  -  V  n  ] 

„j  (3) 

In  formulas  (2)  and  (3) ,  the  symbols  are  the  same  as  those 
used  in  common  MHD  books:  Y  is  the  gas  -  heat  ratio,  N  is  the 
current  sheet's  motion  velocity,  and  n  is  the  current  sheet's 
normal  direction. 

The  kinematic  boundary  conditions  on  the  current  sheet  can 
be  attained  from  the  set  of  equations  in  (2)  by  carrying  out 
integration  on  the  current  sheet: 

[(«  •  v  -  -  0 

n-lH]  -o 

~n(p]  —  p(n  •  V  —  AW V]  +  ±  x  [fi]  X  #}- 0  '  (4) 

[*>(*»  •  V  —  AJ)]  —  o 

In  the  formula,  the  brackets  indicate  the  difference  in  the 
physical  quantity  of  the  two  sides  of  the  current  sheet.  H 
indicates  the  arithmetic  mean  of  H  on  the  two  sides  of  the 
current  sheet . 


We  assume  that 

P  —  Pt*~Mr~’^,y) ,  fi  —  | 

V-V,- («  +  *)*,  Hi- HOC  I  (5) 

a 

In  the  formula,  h=  p  and  /0  indicate  the  pressure  and 

PQ  °  ° 

density  of  the  plasma  on  the  current  sheet,  i  and  K  are  the 
unit  vectors  of  the  x  and  z  directions,  is  the  added 
magnetic  field,  and  H^a  constant.  It  is  very  easy  to  test  and 
verify  that  formula  (5)  satisfies  the  set  of  equations  in  (2) . 
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Taking  formula  (5)  as  the  basic  current,  it  agrees  very  well 
with  the  Wright-Black  solution  [1]  in  the  vicinity  of  the 
current  sheet. 


The  magnetic  field  on  the  left  side  of  the  current  sheet  is 


ffx-Hij 


(«) 


In  the  formula,  H^=a  constant,  j  is  the  unit  vector  of  the  y 
direction . 

In  order  for  the  boundary  conditions  of  (4)  to  be  satisfied, 
there  must  be 


8«  *  8.  • 


(7) 


For  convenience  of  discussion,  below  we  will  take  a  coordinate 
system  on  the  current  sheet. 

2.  Equations  and  Boundary  Conditions  of  Small  Disturbances 
We  superposed  a  small  disturbance  on  the  basic  current  of 
the  plasma  and  magnetic  field.  We  assume  that  all  of  the 
disturbance  quantity  are  shaped  like 

}  - 

In  the  formula,  K  and  Kz  are  wave  numbers.  We  can  then  finally 
obtain  an  equation  which  satisfies  the  disturbance  quantity  v^ 
in  the  plasma: 
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fv.  _  | _ rpUai1 _ 

Jx*  +  rpvKV) 

4* 

+  •  ,*«  •  h  )/{(^+-^^)  [(p*»o +(*;+*».) 

x(r/w^1+-^^«(p#a,J+r^;))]}}  ^ 

~{*S  +  K\  + - - 

rpopo"1  +  e**(A^l«,  +  rPoKl) 


+ - ^ ^ - (*apD«»1Ki(K}  +  Ki)  #*• 

rpopgto1  +  **'(po to1  +  rpoKl) 

4* 

•  a)/{(p»^  +  [(^y  +  (*}  +  «)  (r**w 

+  +  r*KJ))] }}  -  0  (8) 


After  finding  “v  ,  “v  and  v  can  be  found  from  the  following 
x  y  z 

formulas : 


rpt***7  +  ^'(pom1  +  rA#’.)  »'*»  ^ — plauPtKfZ, 

(p^)»  +  (KJ  +  JP.)  {rp«w^  -*-  '*'(w»J+r*Ki)} 

_  r*  {p^  +  X,  %  -  P*iK.  *»(Kl  +  *?)}*- 

(p^)’  +  (KJ  +  «)  {rp+*J  +  **'(p**’  +  r**’.)} 


(9) 


(10) 
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Each  of  the  other  quantities  can  also  be  found  from  the 
disturbance  equation. 

Now,  we  will  find  a  solution  for  the  disturbance  magnetic 
field  on  the  left  side  of  the  current  sheet. 

The  disturbance  magnetic  field  has  potential  and  therefore 

(U) 

The  flux  conservation  of  the  disturbance  magnetic  field  is 
given  as 


W-0  (12) 

can  be  written  as 

<P  -  *  «*'*  •  e~+‘ *”+**.',  K\  -  K\  +  K\,  (13) 

When  substituted  into  formula  (11),  we  can  find  H^. 

For  the  linearization  of  the  kinetic  boundary  conditions  of 
(4)  and  the  kinematic  boundary  conditions  of  (3)  needed  to  fulfill 
the  current  sheet  x=0,  after  eliminating  each  disturbance  quantity, 
we  could  obtain  the  v"  conditions  needed  to  fulfill  x=0: 

X 

[iW  +■  p-  K\j  jr/vsw2  +■  iPV-  (p^ta2  +  +  w  {[iW 

% 

+  kij  rptK\  +  k'J  (pgw2  +  rpoKDj |  v. | 

"  +  < K i  +  **•)  [  rp* 4-  (w* 

3*  1  I-  4* 

+  (14) 


44 


Another  natural  boundary  condition  is 


When  x  —>  +  00  ,  v^  is  limited  (15) 

3.  The  Magnetic  Field  is  Absent  from  the  Plasma  (H^O) 

When  H^O  is  substituted  into  formula  (8)  ,  formula  (8)  is 
transformed  into  a  constant  coefficient  divalent  ordinary 
differential  equation  and  its  solution  is 

Z.  —  <*-**•  A'j  >  0  (16) 


The  the  found  from  the  equation  should  be  satisfactory 

Kl  +  hKt  -  (K’  +  K\ )  -£?*-  X. Hi  (K*  +  K\)  -  0  (17) 

rpt  f  u> 

The  dispersion  equation  is  obtained  from  the  boundary  conditions 
of  (14): 


w*  +  Ki  +  JLH1  k  _ 

r  4.  K. 


PaaKr  —  0 


(18) 


Note  that  the  boundary  conditions  of  (15)  have  already  been 
fulfilled. 

Equation  (18)  shows  that  acceleration  a  causes  the  current 

sheet  to  be  unstable  and  the  greater  a  is,  the  greater  CO  is. 

2  s 

If  a  is  changed  to  a  negative  sign,  then  CO  \  0,  and  the  motion 
is  stable.  Magnetic  field  causes  a  stable  effect. 

If  there  is  only  disturbance  in  the  z  direction,  we  can  find 
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^ I 


(a1  —  a/C, 


(19) 


In  most  situations,  if  K  >  K  ,  we  approximately  take 

2  V 

2 

£0  =aK2  in  formula  (18) ,  and  we  substitute  it  into  formula  (17) 

Kt  -  */k\  +  K\  (20) 

Further,  when  substituted  into  formula  (18)  ,  then 

a*  -  a  -JK\  +  K\  13-ILL  aK\/^K\  +  K\  (21) 

Po  r 

Some  works  [2-4]  have  studied  the  Rayleigh-Taylor  instability 
in  its  experimental  aspect.  Results  have  shown  that  formula  (19) 
qualitatively  agrees  with  the  experiments. 

4.  The  Magnetic  Field  is  Present  in  the  Plasma  (H^O) 

Kruskal  and  Schwarzschild  [5]  have  computed  a  similar 
problem.  They  assumed  that  hx  was  very  small  and  then  transformed 
equation  (10)  into  a  constant  coefficient.  This  cannot  be  done 
in  an  electromagnetically  driven  shock  wave  apparatus. 

If  our  aim  is  only  to  find  the  dispersion  relationship,  then 

from  formula  (14)  we  can  see  that  it  is  only  necessary  to 

dvY  /%/ 

correctly  find  the  value  of  /v  in  x=0. 

Equation  (10)  is  written  as 

+  ^GO  ^  +  BOOS.  -  o  (22) 
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Letting 


1 


then 


+  f  +  ,*(*)/  +  B(«)  -  0 

dx 

When  x  -r>  +  00  ,  formula  (10)  is 


+  K\)Zt  -  0 


(23) 


(24) 


Therefore 

when  x  +  OO  ,  f  (x) 


=  -Y* 


2  2 
y  +  K2 


(25) 


When  x  is  relatively  large,  f  (x)  changes  slowly  and  therefore 
the  conditions  in  (25)  can  be  approximately  changed  to 


when  x=  ,  f  (x)=  /J Jv?"  + 


K 

y  x 


(26) 


Now,  we  use  an  approximation  method  of  equation  (23)  to  solve 
and  fulfill  formula  (26) . 


If 


A(jt)  -  A,  +  AM  B (*)  -  fl.  +  B,(«)  (27) 


in  the  neighborhood  of  x=0,  it  is  assumed  that 
/(.)  -  - K .  +  /,(,)  (28) 

In  formulas  (27)  and  (28),  Kq  is  a  positive  constant,  and  , 
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and  f ^  are  first  order  small  quantities.  After  being  sub¬ 
stituted  into  formula  (23) 


*  Kl  —  AoKt  +  Bt  -  0 


(29) 


and 


+  (<<o  ~  2 AT,)/,  +  S,(*)  -  tfo^.GO  -  0 


(30) 


The  conditions  which  satisfy  f ^  (x)  are 


/,(*)  -  Kt  -  VKl  +  K1.  ^  -  j 


(31) 


After  solving  f ^ (x) ,  we  can  find 


W  +  K\  +  j*(B,0O  -  K„AM)dx |  (32) 

Substituted  into  formula  (14) ,  we  then  obtain  the  dispersion 
relationship. 

If  we  approximately  take  f (x)  |  q=  -Kq,  then  the  dispersion 
relationship  is 


do, 

dx 


/(*)  r  _ 

V, 


-  -K0  4*  * 


f* 
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,  (&)*,  any  k\  ..  -  ,  (r«si  k., 

,  K'  +  r —  ~  poj^o  +  I  (  r  —  l  )p5jio2  x  --2LL 

A«  4»  Kx 

+  (r  —  M  \  f\  rrwo1  +  +  ypaKi)^ 

—  (-  [ rp&ta1  +  (po*oJ  +  rro*-’.)]  {[(wa1)'^;  a|  AT, 

-  ^Kidil  +  KS)  ML'  »}}/{(^).  +  (Kj  k;)  [r(W 


Using  this  formula,  we  can  compute  the  size  of  the  magnetic 
field  required  to  eliminate  Rayleigh-Taylor  instability.  The 
direction  of  the  magnetic  field  is  K  and  this  can  prevent 
z  direction  motion.  Because  of  this,  if  K^=0 ,  the  above  formula 
changes  to 


+  ^-K1,~PtaKt 

4* 


(r-i)(w)a^r  K». 

rpvf**1  +  (pouj1  +  rpoK1.) 

4a» 


Now,  letting  oj  =0,  the  derived  Kq  is 

K\  +  -  ILL  ifssl  -  wK%  -  0  (35] 

4*  T  pi 

2 

In  formula  (29),  letting  u  =0,  the  derived  Kq  is 


*•  — k.  + '3"T^"aT  (Hy*4' 


After  being  substituted  into  C35)  and  simplified,  we  then  have 


2  x  L~ 1  ( -L)1  ■  friJr  -  J_  -  0 

8«*  V  2r  \*./  («})»  ff. 


Fig.  3  is  the  curve  of 

a  disturbance  in  any  z 
is  only  necessary  for 


1_ 

P, 


K)2 

~87f 


and 


h 

K 


Results  show  that  for 


O  '  z 

direction  to  be  able  to  become  stable  it 


i  c/sy  >  r 
p*  8*  r  —  l 


(38) 


Fig.  3 
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MEASUREMENT  OF  ELECTRON  DENSITY  BEHIND  A  STRONG  SHOCK  WAVE 
USING  3  cm  MICROWAVE  TRANSMISSION 

by  Zhu  Naiyi  and  Li  Xuefen 

(Institute  of  Mechanics,  Academia  Sinica) 

Abstract 

This  article  presents  the  working  principle  of  a  3  cm 
microwave  transmission  method.  It  was  used  to  measure  the 
electron  density  behind  strong  shock  waves  in  an  800  mm  diameter 
shock  tube.  The  variation  of  electron  density  with  M  number  is 
obtained,  and  it  is  in  good  agreement  with  other  measurement 
data,  and  is  also  in  good  agreement  with  theoretical  predictions. 

I.  Testing  Apparatus 

Two  split  glass  windows  are  installed  on  the  test  sections 
of  an  800  mm  diameter  shock  wave  tube  [1] .  The  surface  in  the 
glass  and  the  wall  inside  the  shock  wave  tube  is  flat  so  that 
there  was  no  serious  disturbance  of  the  gas  produced.  A  con¬ 
tinuous  iso-amplitude  microwave  signal  is  transmitted  through 
a  horn  antenna  direction  after  300KC  square  modulation.  From 
the  level  of  the  split  window  opening  it  goes  through  the  shock 
wave  tube  and  is  then  received  by  the  receiving  antenna.  In 
order  to  eliminate  the  strong  electromagnetic  field  interference 
caused  in  the  shock  wave  motion  process,  we  floated  the  entire 
receiving  system  and  used  a  symmetrical  output  circuit  to 
transmit  the  signal  from  the  detector  to  the  oscilloscope. 

II.  Working  Principle 

Microwave  transmission  is  used  to  measure  the  electron 
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density  of  plasma  in  shock  waves  based  on  the  relative  changes 
of  the  microwave  power  behind  the  shock  wave  tube.  The  differences 
between  the  shock  waves  reaching  in  front  and  behind  are  mani¬ 
fested  in  two  major  areas.  The  first  is  before  the  shock  waves 
reach  the  window  opening.  The  microwave  interface  has  weak 
reflection  (similar  in  the  areas  of  the  two  window  openings) 
between  the  surface  inside  the  glass  window  and  the  vacuum.  When 
the  shock  waves  reach  the  opening,  the  microwave  interface  area 
has  relatively  strong  reflection  (when  the  electron  density  of 
the  plasma  is  different  so  is  the  strength  of  the  reflection)  in 
between  the  surface  inside  the  glass  and  the  plasma.  The  second 
area  is  before  the  shock  waves  reach.  When  the  microwaves  go 
through  a  shock  wave  tube  which  has  an  800  mm  long  free  space 
internal  area,  there  is  "background"  power  attenuation.  When 
the  shock  waves  reach,  they  penetrate  through  the  800  mm  plasma 
and  can  produce  an  added  power  attenuation.  Based  mainly  on  this 
added  power  attenuation  we  determined  the  electron  density  value 
in  the  shock  wave. 

1.  Reflection  of  the  Surface  Area  in  the  Window  Opening 

The  reflectivity  of  the  microwave  on  the  1  medium  to  2 
medium  reflection  surface  [2]  is 


0) 


It  is  interrelated  with  the  impedance  of  the  medium  wave  on  the 
two  sides  of  the  reflection  surface 


Under  our  conditions,  the  medium's  magnetic  conduction  rate 
was  usually  U^U0'  For  this  reason 


t\-‘  +•  ii/: 


The  relative  dielectric  constant  £  / £  0  of  the  medium  can  be 
written  as 


(f  /e,)*"  -  «  -  jX 


In  this  formula,  n  is  the  common  refractive  index  and 

is  the  attenuation  index,  c  is  the  velocity  of  light,  co  is  the 

microwave  angular  frequency,  and  a  is  the  attenuation  coefficient. 

The  reflective  index  of  the  two  medium  interface  area 
microwaves,  calculated  on  the  basis  of  power,  is 


,  •  | &  |»  —  "iV  ±  (*i  ~  *i)a 

w  lp,’,1  (»,  +  nt)>  +  (*,  +  Xly 


Because  the  inner  diameter  of  the  shock  wave  tube  is  far  larger 
than  the  wave  length  of  the  microwave,  the  electron  density  in 
the  shock  wave  is  close  to  homogeneous.  Given  that  the  electro¬ 
magnetic  wave  is  transmitted  in  infinitely  large,  homogeneous 
Lorentz  plasma,  we  can  attain  the  n  refractive  index  and 
attenuation  index  as  follows 


u>>  ^ 

l  +  if/ 

'j  «>/> 

y+i 

(  V>‘, 

«3  +  V3/ 

'  2  1' 

to1  +  v3 

W  +  V3 

U>1  > 

|  +  J.  \( 

f  2  ' 

J 

y  +  > 

f  U>i 

»3  +  v3  / 

l+2l' 

\  4“  v*  > 

)  + 1 

W  +  v3 

In  these  formulas,  CO  is  the  characteristic  frequency,  and  v  is 

P 

the  collision  frequency  of  electrons  in  the  plasma 

(8  KT\xn 

—)  '  <N‘  •  +  N.  ■  Q„)  (8) 

and  Q en  are  separately  the  electrons  and  ions,  and  the 

neutral  particle  collision  section.  In  the  working  conditions 
of  this  article,  v=10®,  oo  =10^,  that  is  v.  Generally 
speaking,  the  effects  of  v  can  be  disregarded. 

2.  The  Attenuation  of  the  Microwaves  Penetrating  Through 
the  Plasma 

As  previously  mentioned,  in  the  working  conditions  of  this 
article,  we  can  regard  an  angular  frequency  as  a  plane  electro¬ 
magnetic  wave  propagated  in  an  infinitely  large  medium: 

A  —  A,cxp(jut  —  f  •  L)  (9) 

In  this  formula,  the  transmission  constant  is 


f  -  •  +  if 


a  is  the  attenuation  coefficient,  /&  is  the  phase  constant, 

L =800  mm  is  the  distance  traveled  by  the  microwave.  When  v  <&  Co 

—  r.  [(5)  • •  */•  •  i/‘  -  (5)  j  -  n  [(h)/ -J'-’t]  <»> 


54 


2  2 

In  the  formula,  n  =m  to  /4  'Tfe  is  the  cut-off  densitv  of  the 

ec  c 

microwave,  and  ng  is  the  electron  density  in  the  plasma. 

Based  on  the  total  attenuation  of  the  microwave  power  after 
going  through  plasma,  we  deducted  the  effects  of  interface 
reflection  and  found  attenuation  coefficient  a.  Proceeding  to 
the  next  step,  we  obtained  electron  density  ng. 

III.  Conclusion 

The  working  frequency  of  the  microwave  transmission  employed 

in  this  article  is  f=8.87Gc  and  the  corresponding  cut-off 

12  3 

density  is  nec^lxlO  /cm  . 

In  the  800  mm  diameter  shock  wave  tube  presented  in  fig.  1, 

_2 

when  P^=lxl0  mmHg,  under  several  types  of  M  numbers,  there  is 

typical  recording  of  3cm  microwave  transmission.  When  M  is 

s 

approximately  equal  to  12.5,  absorption  cannot  be  seen;  when 

Mg  is  approximately  equal  to  13.5,  weak  absorption  begins;  when 

Mg  is  approximately  equal  to  16,  abs  rption  is  about  one-half; 

when  M  is  approximately  equal  tc  17.6,  cut-off  begins  to 
s 

appear;  when  M  is  approximately  equal  to  20,  there  is  relatively 
s 

long  time  cut-off. 
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Zdff/cm 
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(1) 

(2) 

m 

IB 

(7) 

(8) 


•  lXlO-'ntmHf,  3cm  aWtSUB* 
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1  m(-  i2.5 
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4  M,  - 
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I  M,  »  20.5fe^«Hiil»it 


Fig. 


_2 

1  When  P1=lxlO  nutiHg, 
f=8.87Gc 


the  3  cm  Transmission  Record 


Key:  1.  Ms=12.5  invisible  absorption 

2.  Ms=13.5  weak  absorption 

3.  Ms=14.2  partial  absorption 

4.  Mg=16  absorption  about  one-half 

5.  Ms=17.6  cut-off  begins 

6.  Ms=18.7  partial  cut-off 

7.  Ms=19.9  relatively  long  time  cut-off 

8.  Ms=20.5  relatively  long  time  cut-off 


SOfm/am 


Fig.  2  When  P^=lxlO~1mmHg,  the  3cm  Transmission  Record 
f=8.87Gc 

Test  measurement  attenuation  coefficient  a  given  in  table  1 
changes  with  the  M  number.  Proceeding  to  the  next  step,  we  can 
attain  electron  density  nQ  behind  the  shock  wave  which  follows 
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the  change  relation  of  the  M  number  (see  fig.  3) . 


M, 

13. 5  j  14.2 

15.9 

it 

17.2 

17.6 

0.  56 

i 

0.50 

1 

0.50 

^3 

Table  1  Attenuation  Coefficient  a  Changes  With  the  M  Number 
(P  =lxlO-^mmHg) 

Key:  1.  a  ( neper /meter ) 


Fig.  3  The  Electron  Density  of  the  Air  Behind  the  Shock  Wave 
Pj=lxlO  ^mmHg 


-  Gilmore  [4],  -  Soviet  Union's  Academy  of 

Sciences,  +  Zhu  Naiyi  [3] ,  O  3cm  transmission, 
A  3cm  interferometer,  #  8mm  transmission 


-2 

Within  the  P.=lxlO  mmHg  given  in  fig.  3,  when  M  =12-22, 

X  s 

this  is  the  peak  value  of  the  air's  electron  density  behind  the 

shock  wave.  The  measurement  results  in  this  article,  which  used 


3cm  microwave  transmission,  are  not  only  identical  with  near 
free  molecule  Langmuir  probe  measurements  but  are  also  in  good 
agreement  with  the  reference  values  presented  for  the  3cm 
microwave  interferometer  and  the  8mm  microwave  transmission. 

The  test  measurement  values  were  relatively  close  to  the 
equilibrium  values  which  were  theoretically  predicted. 

Comrades  Luo  Jun,  Li  Shuqin  and  Li  Lianxiang  participated 
in  part  of  the  testing  work  for  this  article.  Comrades  Ye 
Youzhang  and  Gao  Lingzhi  gave  enthusiastic  aid  in  related  work 
for  this  article  and  here  we  would  like  to  express  our  gratitude 
to  them . 
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MEASUREMENT  OF  ELECTRON  DENSITY  BEHIND  A  STRONG 
SH  I  WAVE  USING  3  cm  MICROWAVE  TRANSMISSION 

Zhu  Naiyi  Li  Xuefen 
{liutitvx*  of  Utckaoico,  Aoodomia  IMm) 

Abstract 

This  article  presents  the  working  principle  of  a  3  cm  microwave  transmission 
method.  It  was  used  to  measure  the  electron  density  behind  strong  shock  waves  in 
an  800  mm  dia.  shock-tube.  The  variation  of  electron  density  with  Ms  is  obtained, 
and  it  is  in  good  agreement  with  the  data  from  Langmuir  probe  measurements,  and 
is  also  in  good  agreement  with  theoretical  predictions. 


1.  The  Non-Linear  Interference  Factor 

Based  on  Nielsen's  [1]  proposal,  the  wing-body  combination 
force  can  be  viewed  as  a  combination  of  each  component's  normal 
force  and  interference  quantity 


{K 


watt 


LK»)J> 


s. 


(0 


In  the  formula,  KW(B)P  and  kb(W)P  are  seParate^Y  the  interference 
factors  of  the  wing  when  there  is  a  body,  and  the  wing's  induced 
lift  of  the  body.  The  lower  sign  P  indicates  the  quantity  of  the 
linear  system.  Kw^Bjp  and  KB(W)P  are  identical  [1] . 

Based  on  the  "leading  edge  suction  analogv"  theory,  the 
normal  lift  of  a  single  wing  is  equal  to  [5] 


Cnw  “  Cnm  sinacoaa  +  (i£VLg  + 

CS. 


*r«-CJr(  1  -■—)/«.* 
%vst  *■  2  (~) 


(c;ry 

* A  ' 


(2) 


In  the  formula,  and  ^vse  are  seParate^y  the  leading  edqe 

vortex  and  side  edge  vortex  lift  derivatives.  If  we  assume  that 
the  "suction  analogy"  theory  is  also  applicable  to  the  wing- 
body  combination,  then  at  this  time  it  is  necessary  to  come  forth 
with  more  non-linear  interference  factors 


"*■  +  tCnwir)C\w »in acoaa  •+■  (£*(»>, » 

Ki<w),ru)£rLMnn2a  +  (XrtD.rst  +  — 


Si 


(3) 


In  this  formula,  the  definitions  of  each  non-linear  interference 
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THE  INTERFERENCE  FACTOR  OF  VORTEX  LIFT  OF  WING-BODY 
COMBINATIONS  AT  SUPERSONIC  SPEED 


by  Yin  Xieyuan 

(University  of  Science  and  Technology  of  China,  Department 
of  Mechanics) 

During  the  1950's,  Nielsen  and  others  proposed  the  concept 
of  "interference  factors"  (reference  [1])  to  compute  the 
aerodynamic  interference  among  wing-body  combination  components. 
Their  method  has  become  the  basis  for  presently  commonly  used 
engineering  computation  methods.  During  the  1960's,  comrade 
Ji  Chuchun  and  others  synthesized  and  developed  the  work  of 
Nielsen  [1]  and  Lebedev  [2]  and  brought  forth  a  set  of  computa¬ 
tion  handbooks  [3]  which  have  expanded  the  knowledqe  of  and 
been  widely  used  by  aerodynamic  design  workers.  Yet,  what 
Nielsen  presented  was  linearized  lift  "interference  factors" 
which  are  suitable  for  the  small  attack  angle  range.  At  presenr , 
there  is  still  no  ready  method  which  can  be  followed  for  the 
medium  to  large  attack  angle  range  because  there  is  interference 
produced  within  the  non-linear  lift  by  body  vortex  and  wing 
vortex.  Although,  very  recently,  Nielsen  [4]  brought  forth  the 
concept  of  the  non-linear  interference  factor,  yet  he  still 
used  a  numerical  solution  and  this  does  not  give  a  concrete 
computation  method  for  the  interference  factor.  Other  authors 
have  brought  forth  computation  methods  for  wing-body  combination 
attack  angles  but  most  also  use  numerical  methods.  Based  on  the 
"leading  edge  suction  analogy"  and  "upwash  theory",  this 
article  brings  forth  a  computation  theory  for  the  non-linear 
interference  factor  and  obtains  corresponding  formulas  and 
curves.  It  can  expand  the  application  range  of  the  current 
handbooks  and  is  convenient  for  use  in  design  departments. 
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1.  The  Non-Linear  Interference  Factor 

Based  on  Nielsen's  [1]  proposal,  the  wing-body  combination 
force  can  be  viewed  as  a  combination  of  each  component’s  normal 
force  and  interference  quantitv 


s>  “  +  {Kwni*  +■ 


f } Cnw° 


S, 


(») 


In  the  formula,  KW(B)p  and  Kb(W)P  are  seParately  the  interference 
factors  of  the  wing  when  there  is  a  body,  and  the  wing's  induced 
lift  of  the  body.  The  lower  sign  P  indicates  the  quantity  of  the 
linear  system.  Kw(B)p  and  KB(w)p  are  identical  [1]. 

Based  on  the  "leading  edge  suction  analogv"  theory,  the 
normal  lift  of  a  single  wing  is  equal  to  [5] 
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In  the  formula,  K,,r  „  and  K  are  separately  the  leading  edge 

V  L»£,  Voil 

vortex  and  side  edge  vortex  lift  derivatives.  If  we  assume  that 
the  "suction  antilogy"  theory  is  also  applicable  to  the  wing- 
body  combination,  then  at  this  time  it  is  necessary  to  come  forth 
with  more  non-linear  interference  factors 
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In  this  formula,  the  definitions  of  each  non-linear  interference 
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factor  are  separately 
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The  lower  signs  VLE  and  VSE  separately  indicate  the  quantity 
caused  by  the  wing's  leading  edge  vortex  and  side  edge  vortex. 
If  the  normal  force  distribution  c(t)  )cnv(r)  )  produced  along 
the  spanwise  direction  vortex  is  known,  then 


CflW(BUBLM  “  j]  (4) 


rj  =  X-  and  —  are  separately  the  spanwise  and  chordwise 

Sm  Cr 

dimensionless  distances, 


UCr 


are  separately  the  dimensionless  radius  chord  length  and  wing 
area  (below  we  will  omit  all  of  the  signs).  See  fig.  1. 
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Fig.  1  Diagram  of  Geometric  Measurements 

Reference  [5]  pointed  out  that  the  long  and  thin  triangular  wing 
can  approximately  take  the  vortex  produced  normal  force  as  a 
linear  distribution.  Because  this  article  is  mainly  concerned 
with  the  wing-body  combination  interference  factors,  if  we  take 
the  leading  edge  vortex  produced  normal  force  for  most  of  the 
wing  surfaces  as  having  linear  distribution,  then 

e(.v)ea*(.v)  ”  *(*1  —  Osin1*,/ cos  A,  (5) 

(rj-r)  is  the  spanwise  distance  starting  from  the  outer  exposure 
of  the  wing  root  chord,  a^  is  the  local  attack  angle  because  at 
this  time  the  wing  is  located  in  the  upwash  flow  field  caused  by 
the  body.  If  we  consider  that  the  length  of  the  body  is  sufficient 
and  overlook  the  nose  effect,  then  the  local  attack  angle  dis¬ 
tribution  on  the  symmetrically  level  surface  in  the  cross  flow 
plane  is 
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(6) 


After  substituting  formulas  (5)  and  (6)  into  formula  (4) ,  we 
can  obtain 


Km,hVLt  -  cTTTy  {jr'--Jr2  +  2r-7 
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The  kinematic  conditions  on  the  body  surface  should  be 
sufficient,  that  is,  the  surface  should  be  a  stream  line.  For 
this  reason,  we  should  place  a  corresponding  vortex  image 
system  on  the  body.  Therefore,  the  added  normal  power  of  the 
vortex  produced  normal  force  derived  on  the  body  is 


Because  d  r)  ^  -  ^ydr)  /  c  (  T  i>  cnv  l  1  ^ =c  ( 1  )  cnv  (  1  )  .  therefore 


After  operation  we  obtain 
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Based  on  the  method  found  in  reference  [5] ,  the  side  edge 
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outside  partial  £  angle  is  transformed  into  a  part  of  the 
leading  edge  and  when  £  tends  towards  zero,  we  can  compute  the 
side  edge  vortex  produced  normal  force. 
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Therefore 


K with vse  “  (l  +  r*)1  (9) 

“  rJ(l  +  r2)2  (10) 


In  the  same  way,  we  can  find  the  corresponding  center  of  pressure 
position 
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As  for  X^(B)  ( vsE^t  and  XB(W)  ,VSE/ct'  W€  Can  in  the  same  way 
find 


*»(»>. VS*  —  *K»).VS*  _  _1_ 

e,  e,  2 


(U) 


In  the  formula,  x1  is  the  chordwise  distance  beginning  from  the 
side  edge  apex  quantity. 

(7) -(10)  are  then  four  non-linear  interference  factors  and 
(11) -(13)  are  their  corresponding  center  of  pressure  positions. 
Figs.  2-4  show  their  changes  with  (— )  and  from  the  figures  we 
can  see  that  the  changes  of  the  non-Tinear  interference  factors 
tend  to  be  basically  the  same  as  the  linear  interference  factors 

kw(b)p  and  kb(w)p*  when  r  ~*0,  kw(b),vle  and  KB(W),VSE  *  ll 

m 

and  KB(w)VLE  and  kb(w)vse  when  i - 9l'  they  a11  tend 

m 

towards  4.  Fig.  3  is  the  relative  shift  of  the  wing's  leading 
edge  vortex  produced  normal  force  center  of  pressure.  If  the 
vortex  produced  normal  force  center  of  pressure  of  a  single 
wing  can  be  directly  derived,  then  we  can  raise  the  precision  of 

XW(B)VLE/cr' 
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2  The  Non-Linear  Interference  Factors  of  Wing-Body 
Combinations 


-c. 


3  The  Relative  Shift  of  the  Wing's  Leading  Edge  Vortex 
Produced  Normal  Force  Center  of  Pressure  When  There  is 
a  Body 


Fig.  4  The  Wing's  Leading  Edge  Vortex  Produced  Normal  Force 

Center  of  Pressure  — - — — -  and  the  Normal  Force  Center 

x 

of  Pressure  B  ^  < VLE  Qf  the  Leading  Edge  Vortex  Derived 
r 

on  the  Body  When  There  is  a  Body 


2.  Computation  Examples  and  Discussion 

To  test  and  verify  the  accuracy  of  the  results  obtained  in 
this  article,  we  carried  out  calculations  on  the  data  presented 
in  reference  [6]  and  also  made  a  comparison  with  the  calculation 
method  of  reference  [6]  (fig.  5) . 


Fig.  5(a)  Normal  Force  of  Wing-Body  Combination 


0  Rectangular  wing  and  triangualr  wing  tests  [6]  of 


Key:  1.  This  article  (triangular  type) 
2.  This  article  (rectangular  type) 


Fig.  5(b)  Pitching  Moment  of  Wing-Body  Combination 

0  Rectangular  wing  and  triangular  wing  tests  [6]  of 

Key:  1.  This  article  (rectangular  type) 

2.  This  article  (triangular  type) 

Reference  [6]  considers  that  an  equivalent  rectangular  sub¬ 
stitution  can  be  used  for  all  wing  surfaces  and  it  is  only 

necessary  that  are  all  identical.  When  the  surface  area  and 
m 

center  of  surface  position  are  identical,  then  there  are  identical 
normal  force  and  moment.  His  own  test  data  reveals  this  point 
and  the  two  calculated  curves  of  the  rectangular  wing  and 
triangular  wing  given  in  this  article  are  very  close  to  it.  This 
also  confirms  that  the  non-linear  interference  factors  in  this 
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article  can  be  used. 

The  range  applied  in  this  article  is  the  same  as  that  of 
the  "suction  analogy"  theory.  When  the  attack  angle  is  very 
large,  the  leading  edge  vortex  can  "split"  on  the  wing  surface. 

At  this  time,  the  "suction  analogy"  cannot  be  used.  This  is 
especially  the  case  under  transonic  and  supersonic  speeds  when 
the  large  attack  angle  flow  field  is  very  complex.  Reference 
[7]  points  out  that  at  this  time  the  interference  factor  is  the 
derivative  of  the  attack  angle  and  Mach  number.  At  present, 
there  is  still  no  successful  theoretical  method  which  can 
sufficiently  solve  this  problem.  Aside  from  this,  in  this 
article  we  did  not  consider  the  downwash  effect  of  the  nose 
vortex  on  the  wing  surface. 
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Abstract 

BaMd  upon  lending  edge-suction  analogy  and  npwaah  theory  and  following  the 
•oaeept  of  interference  f acton  proposed  by  J.  N.  Nielaon,  a  method  for  predicting 
the  interference  factors  of  vortex-lift  due  to  nsparation  flow  along  leading  (or  aide) 
edgea  of  wing-body  combination  at  subsonic  speed  is  presented. 

Belated  formulae  and  curves  an  obtained.  They  bars  wider  range  of  applicability 
than  those  available  in  current  handbook 
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